beagle<-read.table("beagle_out_20.snps_chr20.bgl.gz.r2",sep="\t")
dim9beagle)
dim(beagle)
impute<-read.table("chr20.one.phased.impute2_fil_chr20_info",sep=" ",head=T)
dim(impute)
colnames(beagle)<-c("rs_id","r2")
impute<-impute[,c(2,3,6)]
impute[1,]
merged<-merge(beagle,impute,by="rs_id")
merged_f<-merged[which(merged$r2 !='NaN'),
)
merged_f<-merged[which(merged$r2 !='NaN'),]
dim(merged_f)
main<-paste("Beagle r2 VS Impute r2,pearson corr=",cor(merged_f$r2,merged_f$certainty,method="pearson"),sep=" ")
merged_f[1,]
summary(merged_f$r2)
summary(merged_f$certainty)
merged_f$r2<-as.numeric(as.character(merged_f$r2))
summary(merged_f$r2)
main<-paste("Beagle r2 VS Impute r2,pearson corr=",cor(merged_f$r2,merged_f$certainty,method="pearson"),sep=" ")
main
pdf("beagle_r2_vs_impute_r2_chr20.pdf")
plot(merged_f$r2,merged_f$certainty,xlab = c("beagle r2"),ylab=c("impute r2"),main=main,pch=".")
dev.off()
pdf("beagle_R2VSimpute_R2.pdf")
plot(merged_f$position,merged_f$r2,ylim=c(0,1),xlab = c("POSITION"),ylab = c("R2 value"),main=c("beagle r2 (red) vs impute r2 (green)") ,pch=".",col=c("red"))
points(merged_f$position,merged_f$certainty,pch=".",col=c("green"))
dev.off()
multhist <- function(x,main,x.loc,...) {
  new.x <- lapply(x,hist,plot=F)
  tmp <- lapply(new.x,function(x) x$breaks)
  max.break <- lapply(tmp,function(x) max(x))
  max.break <- which(unlist(max.break)==max(unlist(max.break)))[1]
  breaks <- new.x[[max.break]]$breaks
  new.x <- lapply(x,hist,plot=F,breaks=breaks)
  y.max <- max(unlist(lapply(new.x,function(x) max(x$counts))))
  barplot(do.call(rbind,lapply(new.x,function(x) x$counts)),beside=T,
          space=c(0,0),names.arg=signif(new.x[[1]]$mids,2),
          axes=T,legend.text=T,main=main,args.legend=list(x=x.loc,y=y.max),...)
}
pdf("Multiple_hist_chr20.pdf")
multhist(x=list(Beagle_r2= merged_f[,2],Impute_r2= merged_f[,4]),main="Beagle_r2 VS Impute_r2", x.loc=20)
dev.off()
savehistory("plot_chr20.R")
